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COMPUTER-AIDED  ENGINEERING  OF 
SEMICONDUCTOR  INTEGRATED  CIRCUITS 


ABSTRACT 

Economic  procurement  of  small  quantities  of  high  performance  custom 
integrated  circuits  for  military  systems  is  severely  impeded  by  inadequate 

I process,  device  and  circuit  models  to  facilitate  accurate  computer-aided 

j 

design  at  low  cost.  The  salient  objective  of  this  program  is  to  formulate 
physical  models  of  fabrication  processes,  devices  and  circuits  to  allow 
total  computer-aided  design  of  custom  large-scale  integrated  circuit  sub- 

I systems  to  reduce  development  cycle  time  and  cost.  The  basic  areas  under 

investigation  are:  (1)  ion  implantation  and  diffusion  of  dopants,  (2)  thermal 
oxidation,  (3)  chemical  vapor  deposition  of  silicon,  and  (4)  device  simulation 
and  statistical  circuit  modeling.  A brief  discussion  of  the  results  of  the 
first  six  months  of  the  second  year  of  the  program  is  presented  here. 

f 

I 


INTRODUCTION 


In  military,  industrial,  comnercial,  and  consumer  applications,  fre- 
quently there  is  a great  necessity  for  "customizing"  the  design  of  an 
integrated  circuit  to  fulfill  the  critical  needs  of  a specific  system  or 
class  of  systems.  A major  barrier  which  prevents  the  economic  production 
of  small  quantities  of  high  performance  custom  integrated  circuits  is  the 
cost  of  design.  That  is,  the  initial  cost  of  designing  optimum  fabrication 
processes,  device  structures,  and  circuit  configurations  is  prohibitively 
large  because  of  the  amount  of  empirical  human  effort  which  must  be  in- 
vested. The  root  cause  of  this  problem  is  a glaring  lack  of  adequate  process, 
device,  and  circuit  models  and  accompanying  computer-aided  design  techniques 
to  ease  the  burden  of  custom  design.  The  salient  objective  of  this  program 
is  the  development  of  new  basic  models  for  integrated  circuit  processes, 
devices  and  circuits,  which  will  permit  accurate  prediction  of  the  character- 
istics of  a monolithic  structure  on  the  basis  of  its  proposed  process  param- 
eters. These  models  will  serve  as  the  basis  for  economic  computer-aided 
design  of  custom  integrated  circuits. 

The  four  key  generic  integrated  circuit  fabrication  processes  which 
are  being  investigated  are  (1)  ion  implantation  and  diffusion,  (2)  thermal 
oxidation,  (3)  chemical  vapor  deposition  of  silicon,  and  (4)  device  simulation 
and  statistical  circuit  modeling.  This  report  describes  the  progress  which 
has  been  made  in  the  first  six  months  of  the  second  year  of  this  program  in 
all  four  areas. 
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1.  ION  IMPLANTATION 
J.F.  Gibbons,  A.  Chu,  K.C.  Saraswat 


1 .1  Objective 

Our  general  objective  is  to  produce  process  models  for  the  following 
aspects  of  ion  implantation  in  Si  that  are  necessary  for  custom  IC  fabri- 
cation and  which  are  currently  empirical  rather  than  analytical. 

(1 ) Range  and  distribution  calculations  in  multilayer  targets 
based  on  the  Boltzmann  transport  equation. 

- implanted  ion  profiles 

- damage  distribution 

- recoil  implants  ("knock-ons") 

(2)  Diffusive  annealing  behavior  of  implanted  dopants  (multi - 
stream  model ) . 


- profiles  resulting  from  interaction  of  dopants/vacancies/ 
dopant-vacancy  complexes 


- radiation-enhanced  diffusion 

(3)  Study  of  damage. 

- traps  produced  by  knock-ons 

- residual  damage  after  annealing 


1.2  Summary  of  Principal  Accomplishments  to  Date 

During  the  past  two  years  (FY  76  and  77)  the  following  important  contri- 
butions have  been  made. 

1.2.1  Development  of  the  Boltzmann  transport  equation  approach  to 
analytically  predict  range,  damage  and  knock-on  distribution  of  heavy  ions 
(such  as  arsenic),  in  Si,  insulators  (Si02,  Si^N^)  and  a multilayer  combination 
of  them.  Using  the  widely  used  LSS  approach  it  is  not  possible  to  do  all 
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the  above  calculations. 

1.2.2  A multistream  diffusion  model  has  been  developed,  which  apart 
from  simple  thermal  diffusion  of  dopants,  can  also  predict  diffusive  annealing 
behavior  of  ion  implanted  dopants  and  proton  enhanced  diffusion.  The  case 
of  boron  has  been  investigated  by  taking  into  account  boron  ions,  positively 
charged  vacancies  and  neutral  boron  vacancy  pairs  and  their  generation  -- 
recombination  kinetics  to  predict  the  resulting  impurity  profile  of  boron 
implanted  in  the  dose  range  of  10^^  to  10^^  ions/cm^  and  annealed  at  temper- 
atures of  800°C,  900°C  and  1000°C.  Transmission  electron  microscopy  studies 
were  conducted  to  relate  the  diffusion  mechanism  of  boron  to  the  negatively 
charged  vacancies  and  boron  vacancy  complexes.  Using  these  studies  the 
phenomenon  of  precipitation  of  boron  at  higher  dose  implants  was  investigated. 
The  generalized  diffusion  model  can  predict  electrically  active,  inactive 
and  precipitated  boron  as  a function  of  time  and  temperature  of  diffusion. 

A discussion  of  results  obtained  during  the  past  six  months  is  given 
in  Section  1.3. 

1 .3  Detailed  Discussion  of  Results 

1.3.1  Calculation  of  Annealed  Impurity  Profiles  of  Ion  Implanted 
Boron  into  Silicon 

A multi-stream  diffusion  model  is  proposed  for  the  calculation 
of  the  annealing  behavior  of  boron  that  is  ion  implanted  into  silicon  at 
room  temperature  and  subsequently  annealed  [1.1, 1.2].  This  model  is  capable 
of  predicting  both  the  redistribution  and  the  electrical  activation  of  boron 
during  the  anneal,  as  a realistic  model  should.  The  calculated  results 
compare  very  well  with  extensive  experimental  data  reported  in  the  literature. 
The  comparison  includes  samples  that  are  implanted  at  room  temperature  with 
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boron  in  the  dose  range  from  10^^  to  10^^  ions/cm^  and  subsequently  annealed 
in  the  temperature  range  from  800°C  to  1000°C.  This  range  of  dose  and 
annealing  conditions  includes  both  the  typical  applications  of  ion  implantation 
as  it  is  applied  in  the  fabrication  of  devices  and  the  unconventional  cases 
of  high  dose  implants  and  low  temperature  annealing. 

The  annealing  model  is  in  essence  a diffusion  model  extended  to 
include  the  following  situations: 

(1)  precipitation  of  boron,  when  the  boron  concentration  level 
exceeds  the  solid  solubility  limit,  in  high  dose  implantation 
cases 

(2)  trapping  of  boron  by  the  strain  fields  associated  with  dis- 
location dipoles  that  form  during  annealing  at  low  temperature 

In  both  cases,  the  presence  of  this  irmobile  and  electrically  inactive  boron 

•» 

alters  significantly  the  redistribution  and  electrical  activation  behavior 
of  boron  during  the  anneal. 

The  presence  of  boron  precipitates  makes  it  necessary  to  include 
four  species  in  the  annealing  model.  These  species  are  substitutional 
(electrically  active)  boron,  boron-vacancy  pairs  (electrically  inactive), 
positively  charged  vacancies,  and  the  aforementioned  immobile  boron.  The 
electrically  active  boron  is  assumed  to  diffuse  substi tutionally  by  means 
of  random  encounters  with  neutral  vacancies.  The  boron- vacancy  pair  is 
assumed  to  diffuse  much  more  rapidly.  The  principal  interaction  among  the 
species  is  the  reaction  of  active  boron  with  positive  vacancies  to  form 
BV-pairs. 

The  parameters  in  the  model  are  a set  of  diffusivities  and 
lifetimes  for  various  species  and  interaction.  These  parameters  are 
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estimated  from  the  examination  of  each  diffusion  and  each  interaction.  The 
values  of  the  most  important  parameters  obey  simple  activation  energy 
relations.  In  other  cases,  the  modeling  was  either  not  carried  to  the  extent 
necessary  as  to  make  the  temperature  dependence  apparent,  or  else  the  param- 
eters represent  the  composite  effects  of  very  complex  interactions.  In  all 
cases,  the  mathematical  model  and  the  parameter  set  are  universal  in  the 
sense  that  the  same  equations  and  parameters  can  be  used  for  the  prediction 
of  ordinary  diffusion,  proton  enhanced  diffusion,  and  the  annealing  behavior 
of  ion  implanted  boron  at  room  temperature  under  conditions  that  include 
precipitation  effects  and  low  annealing  temperature  anomalies.  Furthermore, 
in  most  cases  of  practical  importance,  high  dose  and  low  annealing  temperature 
conditions  are  of  limited  interest  and,  as  a result,  the  three  stream  dif- 
fusion model  with  an  appropriate  subset  of  parameters  is  sufficient.  This 
reduced  model  is  still  very  important  and  attractive  because  once  the  param- 
eters associated  with  high  dose  and/or  low  temperature  anomalies  are  deleted, 
the  remaining  parameters  have  simple  activation  energies  and  are  determined 
by  the  specification  of  a single  variable  — the  temperature. 

1.3.2  Epitaxial  Silicon  Growth  on  Ion  Implanted  Silicon 

Ion  implanted  technology  offers  an  attractive  means  for  fabri- 
cation of  buried  collectors  in  silicon  bipolar  integrated  circuits.  High- 
dose  implants  (10^^  - 10^^  cm'^)  at  fairly  high  energies  (50-300  keV)  are 
required  for  this  purpose,  resulting  in  heavy  damage  in  the  implanted  region. 
Unless  proper  steps  are  taken  it  becomes  difficult  to  grow  an  epitaxial 
layer  free  of  defects  such  as  stacking  faults.  An  investigation  is  presented 
here  on  effects  of  process  variation  on  the  defects  in  the  epitaxial  layers 
grown  from  SiH^  source  on  <100>  and  <11 1>  oriented  wafers  implanted  with 
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high  doses  of  arsenic. 

Boron  doped  silicon  wafers  of  5-15  fi-cm  resistivity  were  im- 

15  2 

planted  with  100  keV  As,  at  a dose  of  5 x 10  ions/cm  . During  the  implant 
wafers  were  kept  at  room  temperature  or  liquid  N2  temperature.  The  two 
temperatures  were  chosen  to  distinguish  the  effect  of  deep  disordered 
regions  produced  by  the  high  dose  implants  performed  at  higher  wafer  temper- 
atures [1.3].  At  high  doses  the  implant  layers  become  amorphous.  If  an 
annealing  is  done  at  temperatures  above  500*^0,  the  amorphous  layer  on  the 
single  crystal  substrate  regrows  epitaxially  [1.2].  For  high  anneal  temper- 
atures (900-1000°C)  the  resulting  annealed  layers  are  partly  polycrystalline, 
but  a sequential  anneal  results  in  completely  single  crystal  layers  [1.4]. 
Therefore,  two  anneal  cycles  were  employed  in  our  studies.  (1)  A single 
step,  one  hour  anneal  in  nitrogen  ambient  at  900°C  or  (2)  a sequential 
anneal  --  one  hour  at  400°C  in  followed  by  15  hours  at  550°C  in  N2  followed 
by  one  hour  at  900°C  in  N2-  Following  the  annealing  some  wafers  were  subjected 
to  one  of  the  two  drive-ins  at  1250°C:  (a)  15  minutes  in  O2  followed  by 

165  minutes  in  N2,  resulting  in  0.13  urn  thick  Si02  growth;  (b)  10  minutes 
in  ©2  followed  by  60  minutes  in  wet  O2  followed  by  100  minutes  in  O2  followed 
by  10  minutes  in  N2  resulting  in  about  1 un  thick  Si02  growth.  Control  wafers 
without  any  implant  were  also  subjected  to  drive-ins. 

After  stripping  off  the  oxide,  5 um  thick  epitaxial  layers  were 

grown  in  a horizontal  rf-heated  reacter  at  a temperature  of  1070°C.  SiH^ 

partial  pressure  was  kept  at  1.1  x 10  atm  to  give  a growth  rate  of  0.6 

um/min  and  AsH^  partial  pressure  was  kept  at  2 x 10"^^  atm  to  give  As 
15  -3 

concentration  of  10  cm  in  the  epitaxial  layer.  Prior  to  the  growth  the 
wafers  were  subjected  to  an  HC1  etch  at  1180°C.  On  most  of  the  wafers  0.1 
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un  Si  was  etched  from  the  surface,  however  on  some  wafers  the  etch  thickness 
was  increased  to  0.45  un.  After  the  growth  the  wafers  were  subjected  to 
a one  minute  secco  etch  [1.5]  to  reveal  the  defects  in  the  crystal  structure. 
The  following  conclusions  were  reached: 

(1)  For  identical  processing  the  density  of  defects  (stacking 
faults  and  dislocations)  is  more  for  <111>  orientation  than 
for  <100>  orientation. 

(2)  Sequential  annealing  reduces  the  density  of  defects  compared 
to  one  step  annealing  at  900°C.  This  supports  the  study  [1.4] 
described  earlier  in  this  paper  showing  that  samples  with  single 
step  900°C  anneal  do  not  fully  regain  their  crystal  structure. 

(3)  Wafers  kept  at  liquid  N2  temperature  during  the  implant  show 
much  less  density  of  defects. 

(4)  As  compared  to  the  0.1  ijm  HCl  etch  the  0.45  pm  etch  gives  much 
less  defects,  and  the  crystal  quality  of  the  epitaxial  layer 

is  much  less  dependent  on  the  anneal -drive-in  parameters.  This 
should  be  expected  because  the  thicker  etch  physically  removes 
the  damaged  layer  prior  to  the  growth. 

(5)  Wafers  with  1 pm  Si02  growth  during  drive-in  showed  much  less 
defects  than  with  0.13  pm  Si02  growth.  Here  also  the  difference 
being  the  removal  of  the  damaged  Si  by  converting  it  to  Si02. 
However,  the  advantage  of  this  method  is  that  because  the  dopant 
piles  up  in  Si  during  the  oxide  growth,  it  is  not  lost  unlike 
the  HCl  etch. 

(6)  A combination  of  liquid  N2  substrate  temperature  during  implant, 
sequential  annealing,  1 pm  Si02  growth  during  the  drive-in  and 


L.L1^L  i' 
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0.45  ytn  HCL  etch  gave  the  best,  defect  free  epitaxial  layers 
for  both  orientations.  For  the  0.45  pm  HCl  etch  the  1 pm  Si02 
growth  and  liquid  Ng  steps  were  found  redundant;  only  the 
sequential  annealing  was  found  adequate  to  give  defect  free 
epitaxial  layers  as  long  as  the  drive-in  was  done.  For  the 
shorter  etch  of  0.1  pm  the  liquid  N2  step  was  necessary  to 
obtain  the  defect  free  layers. 

(7)  Without  doing  any  drive-in  it  is  not  possible  to  grow  a defect 
free  layer. 


I 
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2.  THERMAL  OXIDATION 

J.D.  Plummer,  B.E.  Deal,  W.A.  Tiller,  C.P.  Ho,  R.R.  Razouk 

2.1  Objective 

The  overall  goals  of  the  thermal  oxidation  portion  of  this  program  may 
be  stated  as  follows. 

2.1.1  To  achieve  accurate  analytic  prediction  of  oxide  thickness  under 
all  process  conditions  encountered  in  modern  technologies.  These  include 
the  effects  of  high  substrate  doping  levels,  arbitrary  silicon  substrate 
orientation,  the  presence  of  a chlorine  species  during  the  oxidation,  and 
the  use  of  multiple  species  (such  as  O2/H2O  mixtures)  during  oxidation. 

2.1.2  To  achieve  accurate  analytic  prediction  of  oxide  charge  densities, 
also  under  all  process  conditions  important  to  modern  devices.  These  condi- 
tions include  arbitrary  silicon  substrate  orientation,  ambient  conditions 
during  oxidation  and  high  and  low  temperature  anneals. 

2.1.3  To  achieve  accurate  analytic  prediction  of  segregation  and  re- 
distribution effects  at  the  Si-Si02  interface.  This  includes  accurate 
prediction  of  dopant  profiles  in  both  the  Si02  and  Si. 

2.1.4  A further  important  goal  is  the  determination  of  an  atomic  level, 
electrically  accurate  model  of  the  Si-Si02  interface.  Such  a model  is 
believed  to  be  essential  for  achieving  physical  understanding  of  oxidation 
kinetics,  the  origin  of  oxide  charges  and  segregation  and  redistribution. 


This  work  represents  a joint  effort  by  the  Stanford  University  Integrated 
Circuits  Laboratory,  Fairchild  Camera  and  Instriment  Corporation  Research 
and  Development  Laboratory  and  the  Stanford  University  Department  of 
Materials  Science. 
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2.2  Summary  of  Principal  Accomplishments  to  Date 


2.2.1  Complete  characterization  of  <111>  and  <100>  orientation  kinetics 
for  dry  O2  and  wet  O2  between  700  and  1200°C  has  been  done. 

2.2.2  Complete  characterization  of  O2/HCI  oxidation  kinetics  for  0-10% 
HCl  and  T = 900  - 1100°C  has  been  completed. 

2.2.3  Complete  characterization  of  heavily  doped  phosphorus  substrate 
kinetics  for  Np  = 10^^  to  solid  solubility,  T = 800  - llOO^C  and  dry  O2  has 
been  completed. 

2.2.4  An  empirical  relationship  of  the  form 

I = C e’*'^  1 + k C 

^ L^e  I + K^tgpe 

where  Cg^  is  the  chemical  phosphorus  concentration  in  the  substrate  and 
Cq,  kp  and  k2  are  constants  has  been  derived  to  model  the  experimental 
results  of  3. 

2.2.5  Implementation  of  computer  prediction  of  first  order  oxidation 
kinetics  has  been  realized. 

2.2.6  The  kinetics  of  H2O  and  H2O/HCI  mixtures  using  a pyrogenic  system 
have  been  studied. 

The  remainder  of  this  section  describes  the  accomplishments  of  the  past 
six  months. 


2.3  Detailed  Discussion  of  Results 

2.3.1  Thermal  Oxidation  of  Heavily  Phosphorus  Doped  Silicon  -- 


Experimental  Results  and  Empirical  Models 


A brief  review  of  general  oxidation  theory  indicates  qualitatively 
what  may  be  expected  in  the  case  of  heavily  doped  substrates  or  diffused 
regions.  The  macroscopic  model  of  Deal  and  Grove  [2.1]  for  silicon 
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oxidation  by  water  or  dry  oxygen  involves  the  diffusion  of  the  oxidizing  species 
from  the  ambient  through  an  existing  oxide  to  react  with  silicon  at  the 
Si/SiOg  interface.  Analysis  based  upon  requiring  continuity  of  steady- 
state  flux  of  the  oxidant  through  the  ambient/oxide/interface  structure 
yields  the  familiar  linear-parabolic  growth  relationship 


where  the  terms  have  the  same  meaning  as  in  Eq.  2.1. 

For  thick  oxides  and  oxidation  times  long  relative  to  a character- 
2 

istic  time  A /4B,  the  first  term  on  the  left  hand  side  of  Eq.  2.1  dominates, 
resulting  in  parabolic,  diffusion  limited  growth.  Under  these  conditions, 

B is  the  dominant  rate  constant  and  is  given  by 

B = 2Dg^^C*/N^  (2.2) 

where  = effective  diffusion  coefficient  of  the  oxidant  in  the  oxide, 

C*  * equilibrium  concentration  of  the  oxidant  in  the  oxide,  and  = number 

of  oxidant  molecules  incorporated  per  unit  volume  of  oxide  grown.  B exhibits 
an  Arrhenius  temperature  dependence  with  activation  energies  of  1.2  eV  for 
dry  O2  and  0.71  eV  for  H2O,  apparently  reflecting  the  activation  energies 
for  the  diffusivity  of  O2  and  H2O,  respectively,  in  Si02  [2.1]. 

For  thin  oxides  and  oxidation  times  short  relative  to  A^/4B,  the 
second  term  on  the  left  hand  side  of  Eq.  2.1  dominates,  resulting  in  linear 
interface  reaction  limited  growth.  B/A  is  the  governing  rate  constant  under 
these  conditions  and  is  given  by 


1 


i 
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B _ kh  C* 

^ “ k + h 


(2.3) 


where  h = gas  phase  transport  coefficient  for  the  oxidizing  species  from 
the  ambient  to  the  outer  oxide  surface  and  k = Si/Si02  interface  reaction 
rate  constant.  An  activation  of  energy  of  2 eV  is  observed  for  B/A  for 
both  dry  O2  and  H2O  ambients,  apparently  determined  by  the  energy  required 
to  break  a Si-Si  bond  [2. 1,2. 2 ], 

In  an  idealized  interpretation,  then,  factors  likely  to  influence  the 


interface  reaction  rate  should  alter  B/A,  while  those  affecting  oxidant 
diffusion  will  change  B.  Sufficiently  high  impurity  levels  in  the  substrate 
should  be  expected  to  modify  the  interface  reaction.  Similarly,  high  impurity 
content  in  the  oxide  may  affect  diffusivity  of  the  oxidizing  species.  The 
relative  magnitudes  of  such  effects  should  depend  on  the  particular  impurity 
involved  and  its  behavior  during  oxidation. 

Specifically,  phosphorus,  as  with  the  other  commonly  used  donor  impuri- 
ties arsenic  and  antimony,  diffuses  more  slowly  in  the  oxide  than  in  the 
silicon  and  tends  to  segregate  at  the  interface  in  favor  of  greater  phosphorus 
levels  on  the  silicon  side.  As  a result,  a pile-up  of  phosphorus  at  the  inter- 
face to  levels  greater  than  bulk  concentrations  in  the  silicon  may  occur 
during  oxidation,  with  relatively  much  lower  levels  in  the  oxide.  In  addi- 
tion, the  pile-up  should  be  more  substantial  for  H2O  than  for  dry  O2 
oxidation  due  to  even  slower  diffusion  relative  to  oxide  growth  rate  [2.3]. 

Therefore,  oxidation  regimes  dominated  by  the  linear  rate  con- 


1 


stant  should  exhibit  greater  dependence  on  heavy  phosphorus  doping  level  than 
should  regimes  dominated  by  the  parabolic  rate  constant.  Shorter  oxidations 
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which  are  interface  reaction  rate  limited  should  show  greater  variation  with 
heavy  phosphorus  doping  than  should  longer  oxidations  which  are  diffusion 
limited.  The  influence  of  heavy  phosphorus  doping  via  large  B/A  variations 
should  also  diminish  or  saturate  with  increasing  oxidation  temperature,  since 
diffusion-limited  oxidation  dominates  at  shorter  oxidation  times  with  increas- 
ing temperature. 

Overall,  then,  this  qualitative  model  indicates  that  greater  influence 
of  heavy  doping  levels  of  phosphorus  (and  also  arsenic  and  antimony)  should 
be  seen  at  shorter  oxidation  times  and  lower  oxidation  temperatures,  and 
that  these  effects  should  be  greater  for  H2O  than  for  dry  oxygen  oxidations. 

The  validity  of  this  simple  model  for  the  enhanced  oxidation  rates  observed 
over  regions  will  now  be  considered  quantitatively. 

An  extensive  series  of  oxidations  were  done  in  a dry  O2  ambient  at  tem- 
peratures of  800,  900,  1000,  and  1100°C.  The  data  produced  are  shown  in 
Figs.  2-1  and  2-2  as  oxide  thickness  vs.  oxidation  time  curves  for  each  of 
five  heavily  doped  sample  types  and  the  lightly  doped  control.  At  each 
temperature,  the  bottom  curve  represents  the  lightly  doped  control  and  agrees 
well  with  published  data  for  <11 1>  substrates,  while  the  additional  curves 
represent  the  heavily  doped  samples,  in  order  of  increasing  substrate  doping. 
Specific  doping  levels  are  indicated  in  Table  I.  Both  SIMS  and  AES  techniques 
were  used  to  measure  chemical  phosphorus  concentration  spreading  resis- 
tance and  anodic  oxidation  techniques  were  used  to  measure  the  electically 
active  concentrations  Cg^. 

The  data  in  Figs.  2-1  and  2-2  already  lends  qualitative  support 
to  the  predictions  based  on  the  simple  interpretations  of  doping  effects 
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presented  earlier.  Clearly  the  greatest  relative  variations  of  oxide  thickness 
occur  at  the  shorter  times  and  lower  temperatures,  the  regions  in  which  the 
interface  reaction  and  hence  B/A  are  expected  to  dominate,  indicating  that 
the  increasing  substrate  phosphorus  level  has  a very  strong  effect  on  the 
interface  reaction.  At  the  longer  times  and  higher  temperatures,  in  the 
diffusion-limited  regime,  substantially  less  variation  of  thickness  with 
doping  is  seen,  implying  that  oxidant  diffusion  and  the  parabolic  rate  con- 
stant are  influenced  to  a much  smaller  extent  by  the  increasing  substrate 
phosphorus  levels.  Indeed,  a few  selected  oxidations  at  1200°C,  when  the 
characteristic  tine  A /4B  for  lightly  doped  substrate  oxidations  is  only 
about  0.5  minutes  and  therefore  oxidation  is  effectively  entirely  B domi- 
nated, produced  essentially  no  variation  of  oxide  thickness  with  increasing 
phosphorus  doping  that  could  be  resolved  within  the  limits  of  interferometer 
accuracy. 

The  variation  of  the  effective  B/A  and  B values  with  increasing 
initial  substrate  total  chemical  phosphorus  concentration,  representatively 
illustrated  in  Fig.  2-3  for  900°C  dry  O2  oxidation,  exhibits  dramatically  the 
expected  behavior.  B/A,  reflecting  interface  reaction-limited  oxidation, 
increases  sharply  by  more  than  an  order  of  magnitude  as  the  phosphorus  level 
rises  toward  solid  solubility.  On  the  other  hand,  B,  which  represents  the 
diffusion-limited  regime,  is  affected  only  slightly  by  the  increasing  sub- 
strate phosphorus  doping.  Clearly,  the  predictions  of  the  simplified  inter- 
pretations presented  have  been  convincingly  confirmed.  We  now  consider  the 

effects  on  B/A  and  3 in  detail. 

(a)  Linear  Rate  Constant.  The  temperature  dependence  of  B/A 


is  shown  in  Fig.  2-4,  where  log  (B/A)  is  plotted  against  1/T  with  substrate 
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doping  level  as  a parameter.  The  lowest  points  are  the  values  found  for 
the  lightly  doped  control  samples,  while  the  higher  points  correspond, 
respectively,  to  the  increasingly  heavily  phosphorus  doped  substrate  types 
used.  To  a reasonable  approximation,  Arrhenius  temperature  dependence  may 
be  assumed  for  the  effective  B/A  values  corresponding  to  each  substrate 
doping  level,  leading  to  the  least-squares  fitted  lines  and  apparent  acti- 
vation energies  of  the  figure. 

High  impurity  levels  may  produce  the  large  increases  of  effective 
B/A  by  altering  the  activation  energy  and/or  the  pre-exponential  factor.  Appar- 
ently the  basic  limiting  activation  energy  of  approximately  2 eV  associated 
with  the  interface  reaction  is  not  substantially  altered  by  the  high  phos- 
phorus levels,  while  significant  dependence  on  dopant  concentration  may  be 
found  in  the  pre-exponential  factor.  A physical  explanation  for  this  rela- 
tively constant  activation  energy  will  be  given  later  in  this  report. 

Some  indication  of  the  nature  of  this  dependence  may  be  gained  via 
examination  of  the  variation  of  B/A  with  concentration  at  a given  temperature, 

as  was  done  in  Fig.  2-3  for  900^0.  The  sharp  increase  of  the  linear  rate 
constant  suggests  an  exponential  dependence  on  Cg^.,  the  initial  chemical  sub- 
strate phosphorus  concentration,  a possibility  that  is  tested  further  by 
investigating  the  variation  of  log  (B/A)  with  Cgj.  in  Fig.  2-5.  Evidently 
log  (B/A)  eventually  exhibits  roughly  a linear  dependence  on  Co„,  and  hence 

Dt 

B/A  an  exponential  dependence  on  Cg^.,  at  each  given  temperature,  but  only 
after  a stronger  initial  increase.  There  may  be  a slight  temperature  depen- 
dence in  the  magnitude  of  this  initial  increase,  but  to  a large  extent  the 
substrate  background  phosphorus  concentration  dependence  is  very  similar  at 
all  temperatures. 
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Clearly,  many  closed-form  relationships  for  this  Cg^.  dependence  may  be 
used  to  approximate  this  data  to  varying  degrees  of  complexity,  in  the 

absence  of  an  accurate  physical  model  for  that  dependence.  A fairly  simple 

form,  producing  a reasonable  first  order  approximation,  might  be 


^0^ 


(B/A)'  = 1 . 


(2.4) 


-E^/kT 

where  C^e  is  the  published  temperature  dependence  of  the  linear 

rate  constant  for  lightly  doped  substrates  with  E-j  2 eV  [2.1].  Thus,  the 
exponential  eventually  dominates,  but  the  power  law  dependence  produces  the 
stronger  initial  increase.  Optimal  and  ko  values  may  be  found  by 
rewriting  Eq.  2.4  in  the  alternative  form 


log 


(B/A)'  - 1 
n, 


1 


'BC 


= log  k^  + k2  C 


BC 


(2.b) 


Plotting  log 


(B/A)'  - 1 


'BC 


vs.  CgQ  should  therefore  produce  a linear  depen- 


dence with  intercept  log  k.|  and  slope  k2.  for  an  appropriate  value  of  n^ . 

Such  analysis  indicates  that  a value  of  n.|  = 0.5  - 1.0  should  be  used,  with 
= 0.5  giving  particularly  consistent  linear  dependences.  Thus,  the  curves 
plotted  together  with  the  experimental  B/A  values  in  Fig.  2-5  have  been  cal- 
culated using  n^  = 0.5  and  the  resulting  appropriate  k^  and  k2  values.  Agree- 
ment appears  very  reasonable.  The  fitted  k-j  and  k2  values  are  presented  vs. 
1/T  in  Fig.  2-6,  confirming  that  the  temperature  dependence  is  relatively 
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slight  when  compared  with  the  overall  dependence  of  B/A.  In  fact,  approp- 
riate constant  values  of  and  k2  independent  of  T conceivably  could  be 

used  for  all  temperatures.  While  agreement  not  surprisingly  is  not  as  good, 

-11  3/2  -21  3 

using  values  of  k^  = 5.9  x 10  cm  ' and  k2  = 3.0  x 10  cm  gives  reas- 
onable agreement  with  the  experimental  results  over  the  entire  temperature 
range  investigated.  Use  of  Eq.  2.4  with  these  values  for  k] , k2  and  n, 
allows  calculation  of  effective  linear  rate  constants  for  arbitrary  substrate 
phosphorus  doping  levels. 

(b)  Parabolic  Rate  Constant.  B has  been  found  to  be  a much 
weaker  function  of  high  phosphorus  concentrations.  Nevertheless,  it  does 
appear  to  increase  somewhat  with  concentration,  as  seen  in  Fig.  2-7,  where 
log  B is  plotted  vs.  1/T  with  substrate  doping  as  a parameter.  At  each 
temperature,  the  lowest  point  represents  the  lightly  doped  control  while  the 
higher  points  correspond  to  the  heavily  doped  samples.  It  is  apparent  that, 
unlike  the  case  for  B/A,  the  concentration  dependence  of  B varies  signifi- 
cantly with  temperature.  Dopant  level  has  a greater  effect  at  lower  temper- 
atures and  virtually  no  effect  at  1100°C.  Also,  for  a given  doping  level, 
the  net  temperature  dependence  is  not  readily  reconcilable  with  a single 
activation  energy  Arrhenius  relation,  implying  changes  in  activation  energy 
and  a temperature-dependent  pre-exponential  factor,  or  addition  of  a new, 
phosphorus-related  term  with  different  activation  energy.  Fig.  2-8,  present- 
ing log  B vs.  Cgj.  with  temperature  as  parameter,  emphasizes  this  additional 
temperature  dependence.  Further,  the  dopant  concentration  effects  at  a given 
temperature  become  evident.  As  substrate  doping  level  increases  (and  conse- 
quently the  lower  phosphorus  level  in  the  oxide  also  rises),  B first  increases 
rapidly  and  then  begins  to  "saturate,"  increasing  only  very  slowly,  if  at  all. 
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As  with  B/A,  a convenient  closed-form  approximation  of  the  B 
variations  may  be  made.  For  convenience,  no  new  explicit  T dependence  is 
included,  but  the  fitted  "constants"  should  then  reflect  any  dependence  thus 
hidden.  The  initial  rapid  rise  and  subsequent  sluggishness  suggest  a possible 
power  law  dependence  on  dopant  levels.  Initial  substrate  phosphorus  level 
is  used  as  the  variable  for  this  analysis,  though  actual  levels  in  the  oxide 
should  be  more  directly  relevant,  as  discussed  in  [2.4]. 


^0^ 


-E2/kT 


H B 


(2.6) 


-Ep/kT 

where  B^e  is  the  published  temperature  dependence  of  the  parabolic 

rate  constant  for  lightly  doped  substrates  with  E2  = 1.2  eV  for  dry  O2  [2.1]. 

Appropriate  and  n2  values  may  be  extracted  in  a manner  analogous  to 
Eq.  2.C.  This  approach  was  used  with  the  extracted  B values  for  800,  900, 
and  1000°C,  with  surprisingly  good  linear  relations  found.  The  calculated 
curves  in  Fig.  2-8  illustrate  the  very  reasonable  fit  to  the  extracted  B 
values  produced  by  the  above  relation  with  the  k^  and  n2  values  shown  vs. 

1/T  in  Fig.  2-9.  The  expected  additional  T dependence  is  clearly  evident 
in  k^  and  possibly,  to  a lesser  extent,  in  02*  A new  additional  activation 
energy  is  strongly  suggested,  and  this  will  be  explored.  Until  a model 
related  to  basic  physical  phenomena  is  derived  to  explain  these  effects  on 
B,  Eq.  2.6  can  be  used  to  calculate  effective  parabolic  rate  constants  for 
substrates  with  arbitrary  phosphorus  doping  levels. 
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2.3.2  HqO  and  HCl/HpO  Pyrogenic  Steam  Oxidation  Kinetics 

The  kinetics  of  steam  oxidation  of  silicon  over  the  temperature 
range  900°  to  1100°C  was  investigated  using  a pyrogenic  system.  This  system 
involves  the  direct  reaction  of  and  O2  in  the  oxidation  furnace  to  form 
H2O.  It  was  found  that  such  a system  does  not  permit  100  percent  (760  Torr) 
water  to  be  produced.  Therefore,  a ratio  of  H2  to  O2  was  selected  which 
provides  a vapor  pressure  of  water  equal  to  that  obtained  in  a previous 
investigation  where  O2  was  bubbled  through  95°C  H2O  [2.1].  This  H2O  pressure 
is  640  Torr.  ; 

A log-log  plot  of  oxide  thickness  vs.  oxidation  time  is  presented 

15-3  ^ 

in  Fig.  2-10  for  <111>  and  <100>  oriented  n-type  silicon  (Cg  = 10  cm  ).  j 

Data  from  [2.1]  obtained  at  1200°C  for  wet  O2  are  also  included  in  the  figure,  , 

with  little  difference  between  <111>  and  --lOOi-  silicon  being  observed  at 
this  temperature.  Data  were  also  obtained  in  the  pyrogenic  system  at  950° 
and  1050°C  but  are  not  included  in  Fig.  2-10  because  of  space  limitations. 

Rate  constants,  B and  B/A,  were  obtained  and  are  plotted  vs. 

1/T  in  Figs.  2-11  and  2-12.  Activation  energies  are  indicated  and  agree 
reasonably  well  with  previous  results.  In  the  case  of  the  linear  rate 
constant,  B/A,  a constant  ratio  of  1.68  is  obtained  between  the  pre-exponential 
constants  of  <111>  to  <100>  oriented  silicon.  On  the  other  hand,  ihe  B/A 
activation  energies  of  the  two  orientations  (2.05  eV)  are  equal. 

Included  in  Figs.  2-11  and  2-12  are  corresponding  data  for  dry 
O2  oxidations.  The  one  important  observation  is  that  previously  the  <100> 

Si  B/A  activation  energy  data  for  dry  O2  did  not  permit  a straight  line  to 
be  drawn  through  the  scattered  points.  If  the  constant  <111>/<100>  ratio  of 
1.68  from  steam  oxidation  is  assumed,  then  a plot  can  be  forced  through  the 

, J 
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dry  ©2  <100>  data  points.  This  has  been  done  in  Fig.  2-12  and  a reasonable 
fit  for  the  data  results. 

Oxidations  in  H2O  with  additions  to  5 volume  percent  HCl  were 
carried  out  over  the  temperature  range  900°  to  1100°C.  Not  only  did  the 
oxidation  rate  not  increase,  it  actually  decreased.  This  decrease  was  found 
to  be  approximately  5 percent,  the  same  amount  as  the  HCl  added.  It  was 
therefore  concluded  that  the  addition  of  HCl  to  an  H2O  oxidation  ambient 
only  serves  to  dilute  the  oxidizing  species.  On  the  other  hand,  the  presence 
of  the  HCl  provides  a cleaning  action,  and  the  oxides  produced  exhibit  much 
less  contamination  effects. 
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Fla.  2-5.  Linear  rate  constant  vs.  initial  substrate  chemical  phosphorus 
doping  level  with  oxidation  temperature  as  parameter. 
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Fig.  2-8.  Parabolic  rate  constant  vs.  initial  substrate  chemical  phosphorus 
doping  level  with  oxidation  temperature  as  parameter. 
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Fig.  2-12.Dependence  of  linear  rate  constant  B/A  on  temperature 
for  the  thermal  oxidation  of  silicon  in  pyrogenic 
H2O  (640  Torr)  or  dry  O2. 
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3.  CHEMICAL  VAPOR  DEPOSITION  OF  SILICON 
K.C.  Saraswat,  T.I.  Kamins,  R.  Re1f,  M.  Mandurah 


3.1  Objective 

The  overall  objective  of  the  research  in  the  field  of  chemical  vapor 
deposition  of  silicon  is  to  develop  physical  models  in  the  following  areas: 

(1)  Kinetics  of  growth  of  epitaxial  silicon 

(2)  Kinetics  of  dopant  inclusion  in  epitaxial  layers 

- growth  of  uniformly  doped  layers 

- transient  response  study  and  system  transfer  function  of 
the  epitaxial  reactor 

- growth  of  arbitrarily  doped  layers  using  the  transfer  function 
approach 

(3)  Effects  of  existing  layer  dopants 

- outdiffusion  and  autodoping 

- substrate  dopant  effects  on  film  induction  and  defects 

- two  dimensional  effects:  geometry  limitations 

(4)  Polycrystalline  silicon  studies 

- atmospheric  and  low  pressure  deposition  modeling 

- physical  properties  of  the  films 

- electrical  properties  of  the  films 

3 . 2 Summary  of  the  Principal  Accomplishments  to  Date 

2.3.1  An  investigation  has  been  conducted  to  relate  the  steady  state 
gas  phase  dopant  (AsH^.  PH^.  B2Hg  and  Sb(CH2)3)  concentration  to  the  resulting 
doping  density  of  As,  P,  B and  Sb,  in  homogeneous  epitaxial  films  of  silicon 
grown  from  SiH^  and  SiCl^,  in  the  temperature  range  of  1050‘’C  to  IISO^C.  The 
investigation  is  novel  for  AsH^,  B^Hg  and  SbCCH^)^.  For  PH^  the  results 


j 


40 


I 

•i 


agree  with  work  done  by  others. 

3.2.2  An  investigation  has  been  conducted  to  study  the  effect  of  growth 
rate  on  dopant  incorporation  in  the  epitaxial  films.  The  dopants  studied 
were  As  and  Sb  (from  AsH^  and  Sb(CH2)2)-  The  trapping  phenomenon  was  observed 
in  the  case  of  Sb.  Similar  study  on  B and  P has  also  been  started. 

3.2.3  Transient  response  of  the  dopant  (As  and  P)  incorporation  process 

has  been  determined  by  studying  the  doping  profiles  that  resulted  from  step 
changes  in  the  dopant  gas  (AsH^  and  PH^)  flow  during  the  deposition  while 
the  SiH^  flow  was  kept  constant.  The  object  of  this  project  was  to  determine 
relationships  between  the  time  variation  of  the  dopant  gas  flow  and  the 
spatial  variation  of  the  dopant  concentration  in  the  epitaxial  film.  By 
studying  the  step  response  of  the  system  a transfer  function  has  been  calcu- 
lated to  achieve  the  above  objective.  This  approach  should  be  useful  in  i 

fabricating  arbitrarily  varying  dopant  concentrations  needed  for  certain  j 

desired  device  characteristics.  i 

3.2.4  A fundamental  study  is  being  conducted  to  understand  the  following  ' 

aspects  of  the  dopant  incorporation  process. 

f 

(1)  Dopant  gas  transport  from  the  gas  distribution  system  to  the 

surface  of  the  epitaxial  film.  i 

(2)  Dissociation  of  the  dopant  gas  into  different  constituents 
in  thermodynamic  equilibrium. 

(3)  Adsorption  of  these  dopant  constituents  at  the  surface  of  the 

silicon.  j 

(4)  Surface  chemical  reactions.  | 

(5)  Desorption  of  the  products  of  the  chemical  reaction. 

(6)  Dopant  incorporation  in  the  crystal  lattice. 

s 

j 
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Thermodynamic  equations  have  been  written  for  each  of  the  above  six  processes. 
By  solving  these  equations  a generalized  model  is  being  developed  to  explain 
various  aspects  of  the  dopant  incorporation  process  in  epitaxial  films.  With 
the  help  of  experimental  results  obtained  in  1,  2 and  3,  various  parameters 
of  this  model  are  being  calculated.  This  model  explains  both  the  steady  state 
as  well  as  transient  aspects  of  the  dopant  incorporation. 

3.2.5  An  investigation  has  been  started  to  study  the  deposition  of 
silicon  at  low  pressure  in  a hot  wall  reactor.  In  the  initial  studies  con- 
ducted, x-ray  and  TEM  techniques  have  been  used  to  determine  the  physical 
characteristics  of  these  films  and  the  effect  of  heat  treatment  on  their 
properties. 

Section  3.3  summarizes  the  principal  results  of  the  first  six  months  of 
the  second  year  of  this  program. 

3.3  Detailed  Discussion  of  Results 

3.3.1  Transient  and  Steady  State  Response  of  the  Dopant  System  of 
an  Epitaxial  Reactor:  Growth  Rate  Dependence 

A physical  model  for  the  dopant  inclusion  into  silicon  epitaxial 

films  is  being  developed  at  the  Stanford  Integrated  Circuits  Laboratory  [3.1, 

3.2].  The  transient  and  steady  state  response  of  the  dopant  system  were 

studied  at  different  growth  rates  to  gather  information  about  the  fundamentals 

of  dopant  incorporation.  The  experiments  showed  the  presence  of  two  regions 

of  operation  and  provide  insight  leading  to  a better  understanding  of  the 

doping  process. 

The  experiments  were  carried  out  in  a horizontal,  water-cooled. 


2 

rf-heated  epitaxial  reactor  with  an  effective  cross  section  of  27  cm  above 
the  susceptor.  Hydrogen  was  used  as  the  carrier  gas  at  a velocity  of  34 
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cm/sec  (at  room  temperature).  Silane  was  used  as  the  silicon  source,  and 
arsine  diluted  in  hydrogen  was  the  dopant  source.  The  corrected  wafer 
surface  temperature  was  approximately  1070°C. 

The  transient  investigation  involved  the  study  of  doping  profiles 
that  resulted  from  step  changed  in  the  dopant  gas  flow  during  the  deposition 
while  the  silane  flow  was  kept  constant.  The  experiments  reported  here  cor- 
respond to  changes  between  arsine  partial  pressures  which  produce  dopant 

15  15-3 

concentrations  of  roughly  1 x 10  and  3 x 10  cm  . The  transient  experi- 
ment was  carried  out  at  0.14,  0.32,  0.6  and  1.3  ym/min.  The  substrates  were 
phosphorus-doped,  (lOO)-oriented  silicon  wafers  with  resistivities  between 
5 and  10  ohm-cm.  The  resulting  impurity  profiles  in  the  epitaxial  layers 
were  determined  by  capacitance-vol tage  measurements  on  planar  p-n  junctions, 
mesa  p-n  junctions  and  deep  depletion  MOS  structures  [3.3]. 

Typical  results  from  the  transient  experiments  are  shown  in  Fig. 
3-1,  where  the  measured  dopant  concentration  is  plotted  as  a function  of 
distance  from  the  surface  of  the  epitaxial  film  for  a decreasing  step 
change  in  the  dopant  gas  flow  during  the  deposition  [3.1].  The  impurity 
profiles  obtained  for  growth  rates  of  0.6  and  1.3  um/min  are  indistinguish- 
able (solid  curve),  and  both  differ  strongly  from  the  one  obtained  at  0.14 
pm/min,  which  changes  much  more  abruptly  with  distance  (dotted  curve). 

An  excellent  fit  to  the  data  can  be  obtained  by  using  the 
expression  [3.2] 


N(x)  = Nj  + (Np  - Nj)  1 - exp 


for  X > x„ 
— 0 


(3.1) 


where 


N(x)  =•  dopant  concentration  at  a distance  x from  the  substrate 
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Nj  * initial  doping  level  for  x < 

Np  = final  doping  level 

L(g)=  decay  length  (a  function  of  the  growth  rate  g) 

= distance  from  the  substrate  at  which  the  transient  begins 

<-0  * 3T„) 

The  same  expression  also  fits  the  doping  profile  corresponding  to  an  increasing 
step  change  in  the  dopant  gas  flow  [3.2]. 

Figure  3-2a  shows  the  decay  length  as  a function  of  growth  rate. 
Each  decay  length  was  obtained  by  using  data  similar  to  Fig.  3-1  together 
with  Eq.  3.1.  At  0.6  and  1.3  pm/min  the  decay  lengths  are  equal  (see  Fig. 

3-1),  but  at  lower  growth  rates  the  decay  length  decreases  with  decreasing 
growth  rate.  Let  us  call  the  low-growth-rate  region  of  operation  region 
I and  the  high-growth-rate  region  (g  ^0.6  m/min)  region  II. 

The  fact  that  in  region  II  the  decay  length  is  independent  of 
the  growth  rate  is  very  significant.  The  time  that  it  takes  to  grow  a 
thickness  of  epitaxial  silicon  equal  to  one  decay  length  is  given  by  L/g. 

Since  in  region  II  the  decay  length  remains  constant,  the  transient  dies 
out  faster  at  higher  growth  rates.  This,  in  turn,  implies  that  in  this  region 
of  operation,  the  mechanism  responsible  for  the  transient  behavior  is  related 
to  the  growing  surface  and  not  the  gas  phase. 

The  steady  state  experiments  involved  studying  the  effect  of  the 
silicon  deposition  rate  on  the  impurity  concentration  of  uniformly  doped 
epitaxial  layers.  The  arsine  partial  pressure  was  kept  at  6.2  x 10'^®  atm, 
and  the  silane  partial  pressure  was  adjusted  to  give  growth  rates  varying 
from  0.13  pm/min  to  1.2  pm/min  (2.7  x 10  ^ - 2.3  x 10^^  atm).  The  thicknesses 
of  the  epitaxial  layers  were  kept  between  5.5  and  8.1  pm  and  were  measured 
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by  a groove-and-stain  technique.  A four-point  probe  was  used  to  determine 
the  resistivity  of  the  layers.  The  substrates  were  boron-doped,  (100)- 
oriented  silicon  wafers  with  resistivities  ranging  from  1-4  ohm-cm. 

The  results  of  the  steady  state  experiments  are  shown  in  Fig. 
3-2b,  where  the  impurity  concentrations  of  uniformly  doped  epitaxial  layers 
are  plotted  vs.  silicon  deposition  rate.  A simple  model  considering  mass 
transport  through  the  boundary  layer,  physical -chemical  reactions  at  the 
growing  surface,  and  thermodynamic  equilibrium  is  sufficient  to  explain  the 
data.  At  growth  rates  greater  than  approximately  0.6  pm/min  (region  II) 
the  doping  density  has  been  found  here  (Fig.  3-2b)  and  in  other  laboratories 
[34]  to  become  inversely  proportional  to  the  silicon  deposition  rate  as 
indicated  by  a slope  of  -1  on  a log-log  plot.  In  this  region  the  dopant 
inclusion  rate  has  reached  a maximum  and  will  be  shown  to  be  limited  by  me 
rate  of  a physical -chemical  reaction  taking  place  at  the  surface  [3.5].  At 
growth  rates  less  than  approximately  0.6  pm/min  the  doping  density  varies 
less  rapidly  than  g'^  and  at  very  low  growth  rates  (region  I)  tends  to 
become  constant  at  an  impurity  concentration  dictated  by  thermodynamic  con- 
siderations. 


3.3.2  A Model  for  Dopant  Incorporation  into  Silicon  Epitaxial  Films 
A model  has  been  developed  to  describe  the  incorporation  of 
dopant  atoms  into  silicon  epitaxial  films  during  deposition  from  a SiH^-AsH^- 
mixture  in  a horizontal,  atmospheric-pressure,  epitaxial  reactor.  The 
model  explains  both  the  transient  and  steady  state  responses  of  the  dopant 
system  [3. 1,3. 5]  and  can  be  used  to  derive  the  "transfer  function"  of  the 
dopant  system  [3.2],  i.e.  an  expression  relating  the  time  variation  of  the 
gas-phase  dopant  concentration  to  the  spatial  variation  of  the  epitaxial 
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layer  dopant  concentration. 

The  model  considers  processes  occurring  both  in  the  gas  phase 
and  at  the  surface.  The  sequence  of  steps  in  the  gas  phase  is  shown  in  Fig. 
3-3:  (1)  Mass  transport  of  the  dopant  compound  (AsH^)  from  the  reactor- 

tube  entrance  to  the  deposition  region.  (2)  Mass  transport  of  AsH^  from 
the  well -mixed  main  gas  stream  through  the  boundary  layer  to  the  growing 
surface.  (3)  Gas-phase  chemical  reactions,  in  which  AsH^  may  dissociate 
into  several  different  As-containing  compounds,  particularly  ASH2.  The 
sequence  of  steps  occurring  at  the  surface  is  shown  in  Fig.  3-4.  (4)  Ad- 

sorption of  the  As-containing  compound  on  the  growing  surface.  (5)  Chemical 
dissociation  into  As  and  H in  the  adsorbed  layer.  (6)  Surface  diffusion  and 
incorporation  of  adsorbed  As  at  atomic  steps  and  kink  sites  on  the  surface. 

(7)  "Burying"  of  the  incorporated  As  atoms  by  subsequently  arriving  Si  atoms 
during  epitaxial  growth.  (8)  Desorption  of  hydrogen  from  the  surface. 

By  analyzing  the  sequence  of  steps  listed  above,  a model  can  be 
derived  to  describe  mathematically  the  doping  process.  Under  steady  state 
conditions  all  of  the  steps  occur  at  the  same  rate.  Under  this  condition 
the  epitaxial  layer  dopant  distribution  is  uniform.  However,  when  the  dopant 
gas  flow  is  a function  of  time,  the  doping  process  may  be  in  a transient 
state.  During  this  transient,  all  of  the  steps  involved  in  the  doping  process 
no  longer  occur  at  the  same  rate,  and  the  correct  physical  picture  is  properly 
described  by  considering  mass-balance  of  the  As-containing  compounds  at  each 
of  the  steps  in  the  process. 

Considering  conservation  of  As  between  each  of  the  steps  mentioned 
above,  the  following  mass-balance-type  equations  result: 
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1.  Mass  balance  of  AsH^  within  the  deposition  volume: 


net  rate  of  AsH, 
etttering  the 

rate  at  which  Asti.,  leaves 

rate  of  cliange  of  AsH, 
concentration  in  ttie 

itlie  iiwin  gas  stream  towards 

deposition 
volume  {r^) 

1 the  wafer  surface  (r  ) 

' m 

main  (|as  stream  wittiin 
the  deposition  volume_ 

Mass  balance  of  Asll^  above  the  ijas-ceysta  1 interfaci 


rate  at  which 

net  rate  of 

net  rate  of 

rate  of  change  of  gas 

AsH, arrives  at 
the'^gas-crystal 

Asllj  adsorption 

gas  (itiase 

(itiase  Asllj  concentra- 
tion above  the 

on  '^the  growing 

chemical 

interface  (r  ) 

m J 

surface  (r^) 

reaction  (r^) 

gas-crystal  interface 

3.  Mass  balance  of  AsH^  above  the  gas-crystal  interface: 


net  rate  of 

net  rate  of  adsorption 

rate  of  change  of  gas 

formation  of 

of  AsH.,  on  the  surface 

phase  As  11^  concentration 
above  tlK>  gas-crystal 

AsH^  (r^) 

- 

L ° 

interface 

4.  Mass  balance  of 


AsH^  adsorbed  on  the  surface: 


net  rate  of 
adsorption  of 

net  rate  of  chemical 
dissociation  of  AsH^ 

rate  of  change  of  the 
surface  densi ty  of 

ASH3  on  the 

on  the  surface  (r^)"^ 

adsorbed  AsH, 

surface  (r^ 

- 

5.  Mass  balance  of 


Asll^  adsorbed  on  the  surface: 


net  rate  of 

net  rate  of  chemical 

rate  of  change  of  the 

adsorption  of 

dissociation  of  AsH^ 

surface  density  of 

AsH2 

surface  ('"(j). 

on  the  surface  (r^^' 

- 

adsorbed  Asll^ 

6.  Mass  balance  of  As  adsorbed  on  the  surface: 


net  rate  of 

net  rate  of  incorporation 

'•ate  of  change  of  the 

formation  of  As 

of  adsorbed  As  into  step  i 

surface  density  of 

adsorbed  on  tlie 

oi-  kink  sites  on  the 

adsorbed  As 

surface (r^  ^ r^) 

surface  (r^)  J 

— 

7.  Mass  balance  of  As  occupying  incorporation  sites: 


net  rate  of 

rate  of  "burying"  the 

rate  of  change  of  the 

incorporation  of 

incorporated  As  by 

surface  density  of  As 

As  into  step  oi' 

silicon  atoms  (r) 

occupying  incorporation 

kink  sites  on  the 

si  tes 

surface  (r^) 

(3.2) 


(3.3) 


(3.4) 


(3.5) 


(3.6) 


(3.7) 


(3.8) 
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When  mathematical  expressions  are  substituted  in  Eqs.  3. 2-3.8  above,  a 
system  of  seven  first  order  linear  differential  equations  results.  These 
equations  contain  the  mathematical  description  of  the  doping  process  and 
can  be  used  to  explain  (1)  the  steady  state  dopant  concentration  in  the 
epitaxial  film  as  a function  of  gas-phase  dopant  molar  fraction  and  epitaxial 
growth  rate  [3.5],  and  (2)  the  transient  response  of  the  dopant  system  and 
its  growth  rate  dependence  [3. 1,3. 5].  They  may  also  be  used  to  derive  the 
"transfer  function"  of  the  dopant  system  [3.2]. 

One  convenient  way  to  handle  this  system  of  equations  is  to  con- 
struct an  electric  circuit  represented  by  an  analogous  system  of  equations. 

An  R-C  circuit  corresponding  to  the  system  of  equations  derived  for  the 
dopant  system  is  shown  in  Fig.  3.5.  Each  resistor  represents  one  of  the 
mechanisms  occurring  during  the  doping  process;  current  is  analogous  to 
reactant  flow,  and  voltage  is  analogous  to  the  molar  fractions  of  the  dopant 
species.  This  circuit  representation  of  the  dopant  system  provides  insight 
into  the  different  mechanisms  taking  part  in  the  doping  process  and  their 
relative  importance. 
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Fig.  3-5.  Circuit  representation  of  the  dopant  system.  F..icli  resisl.ince 
corresponds  to  one  of  the  steps  in  the  iloping  process. 
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4.  DEVICE  SIMULATION  AND  STATISTICAL  CIRCUIT  MODELING 


4.1  Objective 

The  overall  project  goals  for  this  section  deal  with  device  modeling 
and  statistical  circuit  design.  The  technical  approach  involves  developing 
models  and  tools  for  using  the  models  to  predict  process-dependent  perfor- 
mance and  yield  for  custom  integrated  circuits.  Specific  development  goals 
include: 

(1)  One  and  two-dimensional  device  simulation  codes  which 

- accept  simulated  impurity  profiles 

- predict  statistical  performance 

- reduce  to  usable  circuit  models 

(2)  Analysis  and  extraction  techniques  for  devices  to  determine 

- parameter  statistics  and  correlations 

- parameters  most  appropriate  for  process  control 

(3)  Statistical  circuit  design  tools  which 

- use  area  scaling  information  and  process  variations  as  input 
data 

- incorporate  model  parameter  correlations  directly 

(4)  Design  language  and  overall  simulation  capabilities  which 

- include  process,  device,  circuit  and  wafer  statistical 
variations 

- predict  yields  and  aid  in  specifying  design  centers  based 
on  process  variations 


MUtamUNO  FAOIt  MlaMK 
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4.2  Summary  of  the  Principal  Accomplishments  to  Date 

During  the  past  six  months  accomplishments  on  specific  topics  include 
the  following: 

4.2.1  The  demonstration  of  SUPREM  as  a vehicle  to  accept  models  for 
oxidation,  epitaxy  and  ion  implantation.  The  oxidation  models  discussed 
earlier  for  substrate  doping,  orientation  and  ambient  effects  have  been 
implemented.  A first  order  model  for  ion  implantation  clustering  (precipi- 
tation) is  also  available  as  well  as  coupled  species  models  for  boron  and 
arsenic. 

4.2.2  The  use  of  SUPREM  and  a one-dimensional  bipolar  device  simula- 
tion program  to  predict  bounds  on  the  transport  current  for  observed  profile 
variations  across  the  wafer. 

4.2.3  The  development  of  a statistical  bipolar  transistor  model  which 
predicts  75%  of  all  observed  parameter  variations  using  only  the  parameter 
distribution  for  transport  current.  This  model  requires  a factor  of  ten 
fewer  statistical  variables  and  correctly  predicts  parameter  correlations. 

4.2.4  The  demonstration  of  an  efficient  Poisson  solution  program  for 
two-dimensional  devices  of  arbitrary  electrode  and  doping  profile  geometries. 
Preliminary  results  on  finite  element  analysis  programs  for  time-domain 
carrier  transport  and  bipolar  emitter  sidewall  current  transport  are  dis- 
cussed. 

4.2.5  The  fabrication  and  evaluation  of  test  structures  which  identify 
key  geometry-limiting  effects  for  DMOS,  bulk  CMOS  and  SOS  device  technologies. 

4.2.6  Implementation  of  circuit  simulation  capabilities  for  geometry- 
based  macro  models  and  statistical  process  models.  The  statistical  capabili- 
ties also  include  a facility  for  equation  parsing. 


Further  details  on  each  of  these  topics  is  given  in  section  4.3. 


4.3  Discussion  of  Results 

During  the  past  six  months  the  basic  building  blocks  have  been  demon- 
strated for  constructing  statistical  device  models  and  predicting  first 
order  process  dependencies  and  control  conditions  directly  from  process  vari- 
ables. Results  of  this  work  [4. 1-4.6]  illustrate  feasibility  both  for  bipolar 
and  MOS  technologies  of  the  most  recent  results. 

4.3.1  Process  Simulation 

An  essential  tool  for  efforts  described  below  is  a process  simu- 
lator. This  tool  is  needed  to  merge  technology  specific  processing  informa- 
tion with  the  general  process  modeling  results  described  in  previous  sections. 

A process  simulator  called  SUPREM  [4. 7-4.8]  is  presently  being  developed  to 
meet  these  modeling  needs.  Fig.  4.1a  shows  a typical  program  input  data  and 
output  plot  for  a multi-step  boron  processing  sequence.  Further  details  of 
the  program  structure  and  use  are  given  elsewhere  [4.8].  The  basic  simula- 
tion capabilities  are  shown  in  Fig.  4.1b.  Using  SUPREM,  a multi-step  process 
can  be  simulated  to  obtain  a set  of  impurity  profiles.  Typical  results  are 
shown  in  Figs.  4.2a  and  4.2b  for  a double-diffused  bipolar  process  and  a 
CMOS  field  threshold  control  implant.  The  subsequent  use  of  the  profiles  for 
device  and  statistical  process  design  are  given  below. 

4.3.2  Bipolar  Device  Simulation^ 

A first  generation  capability  for  simulating  bipolar  transistor 
performance  directly  from  process  specifications  has  been  demonstrated  [4.9]. 

As  a result  of  present  efforts  this  work  has  been  extended  to  include  device 

^Device  simulator  in  this  report  means  direct  output  of  device  model  paranjeters. 
Device  analysis,"  as  described  later,  means  two-dimensional  analysis  of  the 
device  physics.  The  output  from  the  device  analysis  programs  must  be  post- 
processed  to  generate  model  parameters. 


simulation  based  on  profiles  from  SUPREM  and  outputs  which  bracket  observed 
device  parameter  variations  [4. 1-4. 2].  Fig.  4.3a  shows  a comparison  of 
nominal  simulated  and  measured  device  parameters.  Fig.  4.3b  shows  a com- 
parison of  typical  parameter  spreads  based  on  measured  process  spread  and  the 
corresponding  spread  obtained  by  simulation.  Fig.  4.3c  shows  the  calculated 
upper  and  lower  bounds  on  the  measured  distribution  for  1^.  The  summary  of 
results  to  date  can  be  stated  briefly.  Using  process  inputs  an  adequate 
prediction  of  nominal  device  parameters  as  well  as  correct  bounds  on  para- 
meter variations  can  be  achieved.  However,  the  device  simulators  presently 
available  are  constrained  to  first-order  formulations.  Hence,  many  of  the 
physical  effects  observed  in  high  performance  devices  cannot  yet  be  predicted. 
For  example,  a sidewall  injection  in  bipolar  devices  and  short-channel  punch- 
through  in  MOS  transistors  are  effects  which  are  not  modeled.  Continuing 
efforts  in  device  analysis,  as  described  in  section  4.3,  are  necessary  to 
predict  parameters  associated  with  two-dimensional  device  efforts. 

4.3.3  Bipolar  Device  Parameter  Correlations 

The  observed  parameter  distributions  such  as  the  one  shown  in 
Fig.  4.3c  describe  only  scalar  relationships.  That  is,  similar  distributions 
for  other  parameters  might  initially  be  considered  as  independent.  However, 
evaluation  of  experimental  data  reveals  the  expected  result  that  device  para- 
meters are  correlated  [4. 3, 4. 4].  Fig.  4.4a  shows  the  plot  of  parameters  1^ 
and  Bp  across  a wafer  [4.4].  The  coefficient  R defines  the  degree  of  corre- 
lation. A perfect  straight  line  fit  gives  a value  of  ±1  indicating  perfect 
correlations;  if  the  data  points  i.e.  within  a circle,  a value  of  0 indicates 
no  correlation.  Fig.  4.4b  gives  a graphical  interpretation  of  R in  terms  of 
the  major  and  minor  axis  dimensions  of  an  ellipse  which  circumscribes  the 
data  points. 
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The  correlation  coefficient  assumes  linear  relationships  of  the  form: 


Y = M X + C 


(4.1) 


For  the  parameters  and  Bp,  first  order  theory  predicts  that: 


^ GUp 


(4.2) 


(4.3) 


where  A is  the  emitter  area  and 


*"^6  " J "D“  ’ 

‘^Base  " 


(4.4a) 


r 

' “IT" 

Emitter  ^ 


(4.4b) 


To  the  extent  that  variations  in  A and  GUp  are  negligible,  it  follows  that 
and  Bp  should  be  perfectly  correlated  based  on  variations  in  GUg. 

In  fact,  the  observed  correlation  coefficient  is  0.55  in  Fig.  4.4a. 
To  understand  the  observed  correlation,  other  experiments  were  pursued.  The 
correlation  plots  for  Ig  on  devices  where  the  area  A v/as  varied,  indicate 
that  the  GUg  variations  account  for  87%  of  the  observed  Ig  variations.  Using 
1/Ig  to  estimate  GUg  and  Bp/Ig  to  estimate  GUp,  it  was  determined  that  the 


correlation  for  GUg  and  GU^  was  .94,  reflecting  the  coupling  between  emitter 
and  base  impurity  distributions.  Although  these  studies  do  not  resolve  com- 
pletely the  physics  of  the  observed  vs.  Gp  correlation,  this  is  not  sur- 
prising. The  first  order  model  neglects  emitter  space  charge  and  surface 
effects.  The  two-dimensional  device  analysis  efforts  will  provide  a more 
complete  basis  for  understanding  the  physics  and  predicting  the  correlations 
directly. 

An  example  of  the  application  of  the  parameter  correlation  studies 
is  illustrated  in  Fig.  4.5a  which  shows  the  correlation  coefficient  matrix 
for  the  normal -active  Gummel-Poon  bipolar  transistor  model  parameters  [4.10] 
assuming  the  linear  form  for  correlations  given  by  eq.  (4.1).  A factor  anal- 
ysis study  of  the  data  has  shown  that  75%  of  all  parameter  variations  can  be 
predicted  by  the  model  shown  in  Fig.  4.5b  [4.3, 4.4].  Plausibility  for  such 
a model  can  be  argued  based  on  the  high  correlation  coefficients  observed  for 
in  Fig.  4.5a.  This  simple  model  has  the  key  feature  that  Ig  controls  all 
other  parameter  variations.  This  simplifies  the  characterization  process  as 
well  as  specifications  for  process  control.  The  dominance  of  1^  in  all  para- 
meter correlations  is  characteristic  of  a phosphorus  emitter  bipolar  tech- 
nology. 

4.3.4  Two-Dimensional  Device  Analysis 

The  bipolar  transistor  example  used  in  the  preceding  sections  is 
based  on  one-dimensional  impurity  profiles  and  device  simulation.  However, 
two-dimensional  analysis  is  essential  for  small  geometry  MOS  and  bipolar 
devices.  Present  efforts  focus  on  three  aspects  of  device  analysis.  First, 
an  efficient  two-dimensional  Poisson  solution  program  has  been  developed  and 
used  to  solve  for  potential  distributions  in  devices  with  a variety  of  elec- 
trode and  doping  density  configurations.  Figure  4.6  shows  the  solution  for 
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the  equipotential  lines  within  a buried  channel  CCD  structure.  The  solution 
method  follows  that  described  elsewhere  [4.11].  Ongoing  work  will  define 
the  transport  limits  for  CCD  structures  based  on  electrostatic  and  dynamic 
charge  effects.  The  second  activity  in  device  analysis  is  aimed  at  two- 
dimensional  carrier  transport  including  time  domain  effects  such  as  charge 
pumping  in  SOS  devices.  The  computer  program  used  for  this  analysis  is  at 
an  early  stage  of  development.  A finite  element  method  is  used  to  conserve 
grid  allocations  [4.12].  Figure  4.7  shows  the  grid  used  for  an  SOS  MOSFET 
example  which  is  the  principle  test  example.  Bench  marks  of  this  finite 
element  analysis  program  for  Poisson  and  transport  solutions  will  be  avail- 
able within  the  next  six  months.  The  third  effort  in  device  analysis  addresses 
minority  charge  storage  in  the  sidewall  of  bipolar  transistors.  This  effort 
has  demonstrated  that  the  sidewall  charge  can  be  estimated  using  triangular 
finite  elements  whose  dimensions  can  be  changed  to  follow  current  flow  paths. 
Error  limits  for  this  quasi  two-dimensional  approach  are  being  studied.  The 
next  objective  of  this  work  will  be  to  directly  incorporate  the  method  into 
a device  simulator. 

4.3.5  Test  Structure  Results 

Two-dimensional  device  analysis  which  predicts  geometry  and 
impurity  dependent  effects  requires  corroborative  experimental  data.  Three 
sets  of  test  structures  for  high  performance  MOS  technologies  have  been  fab- 
ricated and  evaluated.  The  lateral  spread  of  CMOS  channel-stop  diffusions 
are  being  characterized  using  special  FET  structures.  Preliminary  results 
show  that  MOSFETs  with  channel  length  parallel  to  the  channel  stop  diffusion 
and  width  dimension  across  the  diffusion  can  be  used  to  obtain  quantitative 
profile  information  for  diffusions  of  widths  of  3 um  and  less.  A detailed 
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analysis  of  the  data  will  continue  during  the  next  six  months.  SOS  test 
structures  fabricated  by  Hewlett-Packard  have  been  characterized.  Device 
measurements  show  that  for  channel  lengths  less  than  10  pm  both  avalanche 
generation  of  carriers  and  charge-pumping  affect  dynamic  conductance. 

Present  efforts  are  now  focused  on  two-dimensional  device  analysis  to  pro- 
vide the  tools  necessary  to  separate  the  effects  of  the  two  phenomena.  The 
characterization  of  DMOS/VMOS  structures  is  the  final  topic  to  be  discussed 
in  this  section. 

Results  of  a DMOS/VMOS  study  are  a final  illustration  of  process 
modeling  for  high  performance  MOS  technologies  [4.5, 4.6]  . Figure  4.8a  shows 
a cross  section  view  of  devices  which  are  fabricated  simultaneously.  Fig.  4.8b 
shows  vertical  impurity  profiles  for  two  fabrication  runs  which  scale  the  physi 
cal  channel  length  from  a maximum  of  2 ym  to  a minimum  of  .5  ym.  Figure  4.9a 
shows  a plot  of  the  linear  region  gain  constant  vs.  inverse  of  physical 
"channel"  length.  The  deviation  from  a constant  slope  for  channel  lengths 
less  than  1 ym  indicates  the  onset  of  transport  which  violates  the  one-dimen- 
sional and  gradual  channel  approximations.  At  high  drain  voltages  the  behav- 
ior also  violates  conventional  theory  in  that  the  depletion  region  spreads 
primarily  into  the  n-type  epitaxial  region.  Fig.  4.9b  shows  plots  of  observed 
output  conductance  compared  with  conventional  theory  and  calculations  based  on 
a solution  of  Poisson's  equation.  The  plots  include  data  for  devices  with 
three  different  channel  lengths.  Model  parameters  for  theoretical  curves 
shown  in  Fig.  4.9b  and  the  data  for  Fig.  4.9a  are  based  solely  on  spreading 
resistance  data. 

The  structures  and  data  shown  in  Figs.  4.8  and  4.9  represent  pre- 
liminary considerations  for  one  specific  technology.  Channel  length  limita- 
tions and  two-dimensional  carrier  transport  effects  depend  critically  on 
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spreading  resistance  data. 

The  structures  and  data  shown  in  Figs.  4.9  and  4.10  represent  pre- 
liminary considerations  for  one  specific  technology.  Channel  length  limita- 
tions and  two-dimensional  carrier  transport  effects  depend  critically  on 
fabrication  technology  constraints.  Availability  of  structures  to  charac- 
terize technology  parameters  and  small  geometry  effects  are  essential  for 
both  device  analysis  and  statistical  circuit  modeling. 

4.3.6  MacroiTwdel s and  Statistical  Circuit  Modeling 

Circuit  simulation  is  the  primary  synthesis  and  verification 
today  for  IC  design.  However,  with  increased  circuit  complexity,  it  has 
become  increasingly  difficult  for  the  designer  to  efficiently  control  both 
geometry  and  technology  inputs  which  effect  circuit  performance.  For  high 
performance  devices,  it  is  essential  that  both  factors  be  accounted  for 
since  each  can  affect  yield  and  performance.  Two  contributions  are  dis- 
cussed in  this  section  which  enhance  the  designers  ability  to  predict  and 
control  circuit  performance  based  on  direct  input  of  a minimum  set  of  vari- 
ables related  to  device  geoemetry  and  fabrication  process. 

For. purposes  of  multiple  cell  design,  a circuit  simulator  must  pro- 
vide a convenient  format  for  assembling  and  scaling  the  sub-blocks  to  achieve 
accurate  predictions  of  circuit  performance.  Figure  4.10a  shows  a circuit 
topology  where  area  scaling  of  devices  is  essential  to  achieve  optimum  output 
speed  and  voltage  drive  capabilities.  Each  sub-block  in  Fig.  4.10a  contains 
an  area-scaled  pair  of  CMOS  devices  as  shown.  Figure  4.10b  shows  the  circuit 
simulation  input  file  used  to  model  this  circuit.  The  macro-model  allows  the 
area  scaling  to  be  inputted  directly  from  sub-block  specifications.  Internal 
program  computations  will  reduce  the  macro-models  to  tabled  sets  of  para- 
meters for  efficient  simulation. 
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Statistical  modeling  with  direct  coupling  to  process-dependent  para- 
meters is  the  second  which  controls  circuit  performance.  It  was  demonstrated 
in  section  4.3.3  that  frequently  the  device  model  parameters  are  interrelated 
parametrically  based  on  process-dependent  physical  parameters.  Figure  4.11 
shows  a statistical  model  implemented  using  the  notation  described  by  McCalla 
[4.13].  The  data  for  an  RTL  inverter  are  given  as  an  illustration  and  only 
Bp,  6|^  and  are  specified  on  the  MODEL  card.  The  output  can  include  a plot 
of  delay  time  vs.  or  other  circuit  parameters  as  desired.  The  statistical 
cards  shown  in  Fig.  4.11  illustrate  several  features.  First,  only  the  model 

I card  has  changed  from  the  conventional  circuit  simulator  form.  It  is  also 

possible  to  specify  the  distribution  on  separate  cards,  allowing  the  conven- 
tional model  to  remain  unchanged.  Second,  Bp  and  B[^  are  specified  in  terms 
of  EQN  1 which  is  a user-defined  equation.  The  equation  is  defined  in  terms 

I of  two  constants  M and  C which  represent  the  general  form  given  in  Fig.  4.5b: 

I 

parameter  = MI^  + C (4.5) 

/ ^ 

i Third,  each  parameter  equation  references  1^  and  the  value  of  I^  is  taken 

from  the  distribution  PI  according  to  EQN  2.  Hence,  1^  and  its  distribu- 
tions control  all  other  device  parameters  on  the  model  card. 

The  output  can  include  a plot  of  delay  time  vs.  I^  or  other  cir- 
cuit parameters  as  desired.  The  statistical  cards  shown  in  Fig.  4.11  illus- 
trate several  features.  First,  only  the  model  card  has  changed  from  the 
conventional  circuit  simulator  form.  It  is  also  possible  to  specify  the 
distribution  on  separate  cards,  allowing  the  conventional  model  to  remain 
unchanged.  Second,  Bp  and  Bp  are  specified  in  terms  of  EQN  1 which  is  a 
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user-defined  equation.  The  equation  is  defined  in  terms  of  two  constants  f1 
and  C which  represent  the  general  form  given  in  Fig.  4.5b: 

parameter  = MI^  + C (4.5) 

Third,  each  parameter  equation  references  and  the  value  of  is  taken 
from  the  distribution  PI  according  to  EQN  2.  Hence,  1^  and  its  distributions 
control  all  other  device  parameters  on  the  model  card. 

The  results  shown  in  Fig.  4.11  represent  preliminary  results  in 
implementing  a statistical  design  language.  Further  work  will  demonstrate 
feasibility  in  using  this  language  for  analysis  of  a variety  of  statistical 
device  models  and  circuit  design  examples.  These  applications  are  critical 
in  demonstrating  the  overall  capabilities  for  statistical  design.  The  basic 
premise  of  this  effort  states  that  proper  design  centering  can  make  it  possible 
to  achieve  acceptable  circuit  yields.  This  design  centering  requires  accurate 
device  models  which  account  for  parameter  correlations.  In  addition  the 
underlying  process  variations  must  be  understood  so  that  designs  can  be  cen- 
tered within  controllable  limits.  The  proposed  statistical  circuit  simula- 
tion capabilities  will  make  it  possible  to  simulate  spreads  in  circuit  para- 
meters resulting  from  process  variations, 

4,4  Conclusion 

In  conclusion.  Fig,  4.12  shows  in  block  diagram  form  how  the  SUPREM 
program  relates  to  the  research  efforts  to  be  discussed.  SUPREM,  based  on 
fabrication  process  specifications,  generates  doping  profiles  and  associated 
physical  and  electrical  information.  The  objective  of  the  work  discussed  is 
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to  develop  device  and  circuit  simulation  tools  to  use  directly  the  SUPREM- 
generated  data.  Specifically,  device  simulation  and  device  analysis  based 
on  SUPREM  input  can  generate  Poisson  and  transport  solutions  for  one-  and 
eventually  two-dimensional  device  structures.  In  addition,  the  electrical 
device  model  parameters  can  be  generated  for  circuit  simulation.  A key  point 
to  be  emphasized  here  is  that  in  addition  to  determining  model  parameters 
from  the  device  simulation,  the  underlying  process  and  geometry  sensitivities 
can  be  explained  using  the  SUPREM  capabilities.  The  final  block  in  Fig.  4.12 
indicates  that  circuit  simulation  based  on  device  model  equations  is  avail- 
able. This  report  focuses  on  two  aspects  of  circuit  simulation  — geometry 
scaling  and  statistical  circuit  modeling.  These  capabilities  are  essential 
for  understanding  the  distributions  of  circuit  performance  parameters  -- 
propagation  delay  for  example  --  which  result  from  both  process  and 
device  geometry  variations.  The  ability  to  simulate  and  understand  process 
control  variables  has  a direct  relationship  to  yield  and  circuit  reliability. 
The  underlying  assertion  of  the  device  simulation  and  statistical  circuit 
modeling  efforts  is  that  for  custom  IC's  the  availability  of  computer  tools 
will  enhance  reliability,  procurement  speed  and  technology  documentation. 
Reliability  depends  on  a precise  understanding  of  technology  limits  and 
device  sensitivities  with  process  and  geometry  variations.  These  limits  and 
sensitivities  can  be  predicted  using  the  tools  discussed  above. 


65 


pevjro  Ciron't  Modoliiui 

[4.1]  Dutton,  R.W.,  D.A.  Divokar,  A.Ti.  Gotualo.’:  and  S.E.  Hanson,  "Coiro- 
lation  of  fabrication  Process  and  Electrical  Device  Parainolcr 
Variatiotis,"  ILE.F  I EDM  lecli.  Diqcst,  Paper  ?4.6,  Dec.  ID?!),  p|). 
609-613.  ■ ■'  " 

[4.?]  Dutton,  R.W.,  D.A.  Divekar,  A.G.  Gon/alo/,  S.f.  Hansei\,  and  D.A. 
Antoniaclis,  "Correlation  of  fabrication  Process  and  Electrical 
Device  Parameter  Variations,"  to  be  publislied  in  I EEljl^  of  Solid- 
State  Circirits,  Vol . SC-12,  Auej.  1 977.  '*  “ 

[4.3]  Divekar,  D.A.,  R.W.  Dutton  and  W.J.  McCalla,  "Bipolar  IC  Device 
Statistics  - An  Experimental  Studv,"  Proc.  of  IfEf  ISCAS,  Paper 
SIS. 3,  April  19/7,  pp.  56S-S69. 

[4.4]  Divekar,  D.A.,  R.W.  Dutton  and  W.J.  McCalla,  "Experimental  Study  of 
Gumtiel-Poon  Model  Parameter  Correlations  for  Bipolar  Junction 
Tratisistors, " submitted  to  IfEF  J.  of  Solid-State  Circuits. 


[4.5]  Combs,  S.R.,  D.C.  D'Avanzo  and  R.W.  Duiton,  "Characterization  and 
Model  inrj  of  Simultaneously  fabricated  DMOS  and  VMOS  Transistors," 
lEDf-1  Technical  Pi <].c_s_t , 19/6. 

[4.6]  li'Avanzo,  D.C.,  S.R.  Combs,  and  R.W.  Dutton,  "Effects  of  the  Diffused 
Impurity  Profile  on  the  DC  Cfiaracteri sties  of  VMOS  and  DMOS  Devices," 
IEEE  Solid-State*  Circuits,  Auej.  1977. 

[4.7]  Dutton,  R.W.,  A.G.  Gonzalez,  R.D.  Runcj,  and  D.A.  Antoniadis,  "1C 
Process  Etujinccriiu]  Models  and  Applications,"  conduct  or  Si  l^u^n 
1977,  presented  at  the  Spj'in^  Electroclu'inical  Society  Meeting, 
Philadelphia,  May  9-13,  1977,  pp.  910-922. 

[4.8]  Antoniadis,  D.A.,  S.E.  Hansen,  R.W.  Dutton,  A.G.  Gonzalez,  "SDPREM  1 -- 
A Program  for  1C  Process  Modeling  and  Simulation,"  Stanford  Technical 
Report.  SEE  77-006,  May  1977. 

[4.9]  Gonzalez,  A.G.,  S.R.  Combs,  R.W.  Gill,  and  R.W.  Dutton,  "fabrication 
Process  Modeling  Applied  to  IC  NPN  Transistors  Using  a Minicomputer," 
International  Electronics  Convention  1975,  August  25-29,  Sydney, 

Austral  la. 

[4.10]  Nagel,  L.W.  and  D.O.  Pederson,  "Simulation  Program  with  Integrated  Circuit 
Emphasis  (SPICE),"  presented  at  the  19th  Midwest  Symp.  Circuit  Theory, 
Waterloo,  Ont.,  Canada,  April  1973. 

[4.11]  Kennedy,  D.P.,  "Computer  Aided  Engineering  of  Semiconductor  Integrated 
Circuits,"  Research  and  Development  Technical  Report,  ECOM-75-1.344-1  , 

April  1976. 

[4.12]  Barnes,  J.J.,  R.d.  Lomax,  "Finite  Element  Simulation  of  GaAs  MESfETS  with 
Lateral  Doping  Profiles  and  Submicron  Gate,"  IEEE  Trans.  Electron  Devices, 
Vol.  ED-23,  No.  9,  Sept.  1976,  pp.  1042-1048. 


\ 


[4.13]  McCalla,  W.J.,  "Computer-Aided  Circuit  Simulation  Techniques,"  EE  392B 
Course  Notes,  Stanford  University,  Spring  1977. 


ms  FACffi  IS  BSST  QUALITY  FRAClteAlU 
ISOM  OOFY  FUKNilbhliD  TO  ODC  - 


uu 


Figure  4 


Output  Data 


Input  Data 


I'l*  •if^*!***,  MOavtcwt*.  *oo-.**tri 

fic»  M#***!!!*  IINIat*  ■Miavat*. 


INPUT 

FILE 


INPUT 

SCANNER 

ISUPERVlSORl 


(a) 


SUPREM 

MULTI  STEP  PROCESS  SIMULATOR 


PROCESS  MODELS 


>{sURrACE  DEP. 

0 

1 

■ OXYGEN  AMB. 

^ - 

u : 

' NEUTRAL  AMB. 

S 

1 

0 

• EPITAXY  

N 

J H 


r*-» 


3>ilON  IMPLANTATION 


L 


t>j0XIDE  PEPOSITION/CTd7 


IMPURITY 

PROFILES 


=i  bur  PUT  GENERATOR 


OUTPUT 


(b) 


The  SUPREM  program  (a)  input  and  output  examples,  and  (b) 
overall  program  structure. 


the 


IIHIS  PACT  IS  BEST  QUALITY  PRACTICABU 

reOM  COPY  FURNISHliD  TODDC  ^ 


r.llPI<i:H  II  ST  KXt'.lif'l.i;  :il 
sii:f'  II  (>  ciii’siJi.i:  DiT’usninM  ((ll•::iI:r^I(:) 


0.0000  0.00 

O.lli  0.3(1 

0.4li  O.GO 

. - 

DEPTH 

(IJM) 

SIIPRCM 
,,  STEP  « 

TEST  EXflMPl.E  h 1 
4 oxiiiniiow 

(nRY/lJET/DPY) 

~e.3U0ll  0.00  O.IS  0.30  0.4li 

DLI'TM  (DM) 


o.bo  D.vc  o.no 


Figure  4,2.  Two  output  plots  from  SUI'UfM  (n)  .i  boron  redistribution 
ex.iinple,  oiul  (b)  boron  aiul  arsenic  profiles  for  a bi- 
polar process. 


SImulAlrJ 

H'*AAurcd 

l,(M 

• 16 

7.1$  K 10 

7,*',  K 10'“ 

•r 

106.00 

174.00 

*j(n) 

3*.>0 

<7.63 

v^<v) 

20.  JO 

17.67 

c^.,(rr) 

0.391 

0.36 

0.37 

0.31 

0.30 

0.32 

0.721 

0.77 

''y(rsEC) 

231.00 

179.00 

VArlAtion 

f**r 

Stnulatton 

HCARUlcd 

VArlAtion 

XJ» 

4 3.41X 

4 7X 

Spread Ing 

XJnc 

4 j.eex 

4 <.61X 

tenlAtAncA 

•CAturementt 

btltter  ShcAC 
Rculstance 

4 2.93X 

4 3.<X 

on  the  teet 
vAf  er 

flnch^d  Barc 

Sheet  ScAlfttAncc 

419.  KX 

4 79X 

1« 

436. 8X 

4 <0X 

Oevlco 

f 

®or.s 

4 9.2X 

4 8.  Jt 

■caeureeicnta 
on  test  chip 

4 a S6  6 4 7.2 


Figure  4.  3.  Simulated  and  measured  parameters  for  a bipolar  device 
technology:  (a)  nominal  transistor  electrical  param- 
eters, (b)  physical  and  electrical  variations,  and  (c) 
measured  distribution  and  siimilated  bounds  (dashed). 


Figure  4. 


rais  PAffi  IS  BIST  Quttin  psAcneJM! 

jrOM  COf  Y FURlIlSHJiD  TO  B3C 


Pf 


(a) 


CORRELATION  COEFFICIENT 

R.  r(X-X) (Y-Y) 
V^£(X-X)'-  j:(y-y)* 


|R|.|- 


MINOn  AXIS 
MAJOR  AXIS 


(b) 


4.  Electrical  parameter  correlations,  (a)  Scatter  plot 
for  and  0^  with  correlation  coefficient  of  R = 

.55,  and  (b)  physical  interpretation  of  R factor. 


1OT«W«  IS  BEST  QUALITY  nUCnCIBW  LINTAR  COf^RI  t A ( ION 
I®a«OOPYfutu»iw.u-.  ivue^G  COI'ITICIHN  I MATRIX 


Figure  4. 


Xi  = k Xj  ^ c 


ftr 

A„ 

Rn. 

ft 

1= 

f't 

Ilf 

Va 

(■; 

^or  1 

'/•in 

'"in 

A 

10 

.7  3 

76 

55 

15 

-or. 

01 

-.0 

01 

72 

-.15 

.6 

A, 

10 

7C 

.&/ 

.§ 

50 

5? 

57 

• 06 

-M 

- 24 

.25 

.66 

Hpo 

10 

.56 

65 

6 

67 

7 

-.7 

.71 

-.7 

67 

.31 

•i 

to 

.6 

.1 

.41 

?il 

- 51 

.45 

• 55 

41 

.70 

Is 

10 

54 

56 

50 

• 61 

77 

• 41 

55 

50 

Rt 

10 

.06 

-.7 

-30 

(4 

• 70 

60 

04 

"c 

10 

M 

-.75 

.62 

-.72 

.66 

.04 

10 

-42 

.25 

.15 

-.22 

.55 

L_  . 

10 

-46 

23 

-.2 

-.56 

Cf 

1.0 

.45 

.26 

Qocs 

10 

-.7 

-.07 

^Yb 

1.0 

. 21 

mjD 

1.0 

5 

S.9B  Ij  ♦ 1J9.99 

5 

0.19  Ij  ♦ l.I 

Vi 

16.39  Ij  ♦ 23.19 

(K-OHm) 

V 

43.27  1,  + 739.83 

(Ohns) 

30.63  Ij  - 136.85 

**# 

0.06  Ij  ♦ l.I 

0.28  Ij  + 3.86 

(Ofias) 

V 

-10.93  Ij  + 106.53 

(Voltl) 

*T 

13.87  Ij  ♦ 387.67 

OOU) 

^OF.i 

-0.01  Ij  ♦ 0.5933 

(ri) 

Vi 

0.01  Ij  ♦ 0.31 

^Ei 

O.Ol  Ij  ♦ 0.71 

(Voltl) 

(Ij  l«  a aulttple  o(  10 
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Figure  4.7.  - A two-dimensional  finite  element  mesh  for  SOS  MOSFET. 


DMOS 


SHIS  PAGE  IS  BEST  QUALITY  PlUCriOBliS 
JBOM  COPY  FURtilSRED  TO  DDC 


to 

E 

to  10 

S 

p 


o ONE  MICRON  CHANNEL 
o HALF  MICRON  CHANNEL 


o 

2 

O 

o 


.n  .8  1.2 

DEPTH  IN  MICRONS 


Figure  4.8  DMOS/VMOS  test  structure,  (a)  Fabrication  sequence,  and 
(b)  vertical  doping  profile  for  different  channel- 
length  designs. 


/Xi  = /9LE/ZC, 


WIS  FAOB  IS  BIST  QUALITY  PRAmCAto* 

raCttl  COPY  FURNISHED  TO  DDC 


§Le=0.59/x 


CVJ 

> 

< .8 

E 


/Aj-I73v^/ 


.8  1.2 

I/L£  [\/^] 


&Or 


o SIMUlm.w 
— MEASURED 


Le  = I.29  Na=5.3xI0'5 


Lr  = 2.42  Na  = IxIO'® 


Lg'0.6^ 

N^=  5.9x10*' 


MEASURMENT 

o GAUSSIAN  PROFILE 
* o STEP  FUNCTION 


1 25I 


I 


OOGOOOOOOOOO 


5 

Vq  [volts] 


Vq  (Volts) 


Figure  4.9  Electrical  device  characteristics  for  the  VMOS 
structures:  (a)  linear  region  plot  of  B vs  Y, p 
and  (b)  output  conductance  plots  and  compar- 
ison  with  model  results. 


75. 


IBIS  PA(X  IS  BEST  QUALITY  rRACtlfiAtJJL 
raOM  COFI  FURNISHED  TO  DDC  - 


* 


VIN 

1 

0 

PULSE  0 12  lOM  ION  lON  SON  lOON 

XII 

1 

2 

INV  1 "I 

XI2 

2 

3 

TWV  ^ L 

XI.^ 

3 

4 

INV  10  J 1 

CL 

4 

0 

* 


MACRO  1 ? 0 INV  i.O  J 

Ml  8 1 MP  10.'/ 1 10 

M?  2 10  0 KW  10»l“l0  ■=•  1,3  Ot^  lO 

ME  NO  ^ 

* 

MOOEL  MN  NMUS  VTi)=  llf-<  = J0O  TOX=lOOO  ONH=i5Kl5 
MGOEL  MP  HMOS  VT0  = “«;  Lib  = 2b0  rOX=1000  UNfi=lE;ib 

* 

CND 


Figure  4.10.  Macrotnodeling:  (a)  equivalent  circuit  for  CMOS 
inverter  drain,  and  (b)  circuit  simulator  input 
including  device  area-scaling  feature. 


I 


I 


I 


I 


I 


76 


COPY  PUtiAliHfiD  rODiDC 


* 

* simple  ktl  inverter  circuit 

n 

VCC  ? 0 DC  5 

VIN  1 0 PULSE  b 0 2N  i?N  2N  3aN  «0i\ 

* 

KC  2 4 1 
RB  1 3 10 

♦ 

01  430  OMOO 

* 

MODEL  OMOD  NPN  BFM  17b  EQUATION  EQNl  ( 5 . 9P , QMOD : I S , I 39 . 99) 

♦ bRM  2.2  EQUATION  EQNl  ( 0 . 1 9 , Q'-lOD  : I S , 1 . 2 ) 

♦ IS  6.0t-lb  EOUAnON  EQN2  (PI) 


* 

* USE  OF  OU«MY  parameters 


* 

PI  DISTRIBUTION  TABLE  TPl  NOMINAL  MEAN 


* 

EQUATION  EGiN2  (X)  = X * lE-16 


« 

EQUATION  EONl  (M,X,C)  = M * X * 1E16  ♦ C 

* 

TABLE  TPl  CONTINUOUS  DENSITY 

♦ 4.0,0  4.B, 10  b.b, 13  6.4,3 

♦ 7.2,4  B.0,3  d.«,0  9.6,2 

♦ 

END 


Figure  4.  11.  RTL  inverter  circuit  and  statistical  circuit  simulator 
input  format.  EQN  1 defines  parameter  correlation 
equations  and  EQN  2 gives  distribution  data. 


emmOBBBJaSSiam 


77 


IBIS  PAGE  IS  BEST  QUALITY  PRACITCAftLt 
EROM  COPY  FURMSHJEX'  ID  DDC 


Fig.  4.12. 


79 


5.  SUPREM  I - A PROGRAM  FOR  IC  PROCESS  MODELING  AND  SIMULATION 

The  Stanford  University  PRocess  Engineering  Models  (SUPREM)  program  is 
a computer  simulator  capable  of  simulating  most  typical  IC  fabrication  steps. 

The  program  is  designed  so  that  these  steps  can  be  simulated  either  individually 
or  sequentially,  just  as  they  would  occur  during  the  actual  fabrication 
process.  The  output  of  the  program,  available  at  the  end  of  each  step, 
consists  of  the  one-dimensional  profiles  of  all  the  dopants  present  in  the 
silicon  and  silicon-dioxide  materials.  These  profiles  are  displayed  in 
various  formats  including  line-printer  printout,  line-printer  plots,  and 
high-resolution  (Calcomp  type)  plot.  It  is  understood  that,  in  sequential 
step  simulation,  the  output  of  a process  step  constitutes  the  initial  condi- 
tions for  the  following  one.  The  junction  depths  and  sheet  resistances  of  all 
n or  p layers  formed  during  the  process  are  also  calculated. 

The  fabrication-step  simulation  is  based  on  several  process  models. 

The  models  implemented  in  SUPREM  are  the  following: 

(1)  oxidation/drive-in 

(2)  predeposition  through  the  surface  (gaseous  or  solid) 

(3)  epitaxial  growth 

(4)  ion  implantation 

(5)  etching 

(6)  oxide  deposition 

Solid-state  diffusion  is  fully  accounted  for  in  any  of  the  above  models 
involving  high  temperature. 

Communication  between  the  user  and  SUPREM  has  been  made  as  simple  as 
possible.  The  input  file  consists  of  free  format  statements  involving  key 


•fir-  • 

^ x ffilCXDlllO  PiOB  MUJK. 

I 


80 


words  and  numbers.  Typically,  a single  processing  step  can  be  simulated 
with  an  input  specification  involving  less  than  60  alphanumeric  characters. 

Several  physical  parameters  have  been  stored  in  the  program  and 
constitute  the  default  values  used  by  the  models.  They  may  be  overridden 
by  user-specified  parameters. 

The  following  table  lists  the  organizations  currently  using  the  SUPREM 
program. 


PARTIAL  LIST  OF  DISTRIBUTED  COPIES  OF  SUPREM 


Government 


1. 

ECOM,  Fort  Monmouth,  N.J. 

Randy  Reitmeyer 

2. 

National  Bureau  of  Standards 

Tom  Leedy 

3. 

National  Security  Agency 

Bill  Bundy 

4. 

Harry  Diamond  Laboratories 

Brian  Biehl 

DOD 

Contract  Holders 

1. 

TRW  Defense  and  Space  Systems  Group 

John  Choma 

2. 

Hughes  Aircraft  Co.,  Carlsbad  Research 

Center 

Michael  Jack 

3. 

Sandia  Laboratories 

Chuck  Gwyn 

4. 

Northrop  Corporation,  Northrop  Research 

and  Technology  Center 

F.W.  Pfeiffer 

Industry  U.S.A. 

1. 

RCA  Solid  State  Division,  Sommerville,  N.J. 

Chris  Davis 

2. 

Bell  Laboratories 

a.  Allentown,  Penn. 

Paul  Langer 

b.  Murray  Hill,  N.J. 

Dr.  Jacques  Ruch 

c.  Holmdel,  N.J. 

Alexander  Voshchenkov 

3. 

Intel 

Chris  Martensen 

4. 

Prime  Computer  Inc. 

Mike  Payne 

5. 

Texas  Instruments,  Dallas,  Texas 

Dr.  William  Daughton 

6. 

Tektronix,  Beaverton,  Ore. 

Dr.  Doug  Ritchie 

7. 

Signetics,  Santa  Clara,  CA 

Tak  Young 

8. 

Litronix:  Advanced  LSI  Technology  Inc. 

Cupertino,  CA 

Molz  B.  Khambaty 

9, 

Control  Data  Corporation,  Minneapolis,  Minn. 

Chuck  Naber 

10. 

Bell  Northern  Research 

Michael  Caughey 

n. 

National  C.S.S.,  Santa  Clara,  CA 

Mr.  Wong 

12. 

National  Semiconductor,  Los  Altos,  CA 

Darrell  Erb 

13. 

Xerox,  El  Segundo,  CA 

Kerning  W.  Yeh 

14. 

Philips  Laboratories,  Briarcliff  Manor,  N.Y. 

Jerry  S.  Sullivan 

15. 

Solid  State  Scientific,  Montgomeryvil le,  Penn. 

Jeff  Steinwedel 

16. 

Hewlett-Packard,  Corporate  Engineering, 

Palo  Alto 

Hal  Daseking 

a.  Loveland  Instruments  Division 

Dick  Toftness 

b.  Microwave  Semiconductor  Division 

Craig  Snapp 

c.  Solid-State  Laboratory 

Ken  Lisiak 

d.  Integrated  Circuits  Laboratory 

Skip  Rung 

e.  Data  Systems  Division 

Wei  Wu 

f.  Colorado  Springs  Division 

Grant  Somners 
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Table  I (cont'd) 


Industry  Foreign 

1.  Centre  Electronique  Horloger,  Switzerland 

Educational  Institutions 

1.  University  of  Maryland 

2.  University  of  California,  Berkeley 

3.  University  of  Toronto 


Heinz 


Prof. 

David 

Prof. 


W.  Luginbuhl 


H.C.  Lin 
Angst 

C.A.T.  Salama 
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6.  ONE  DAY  SYMPOSIUM  ON  PROCESS  MODELING 

A one  day  symposium  was  held  at  Stanford  University  on  July  8,  1977, 
on  Integrated  Circuit  Process  Engineering  Models.  The  main  objective  was 
to  provide  a sunmary  of  our  recent  advances  in  predicting  and  understanding 
ion  implantation,  diffusion,  oxidation  and  epitaxy.  The  content  of  the 
symposium  reflected  ongoing  research  at  Stanford  under  DARPA  sponsorship. 

In  all,  six  presentations  were  given  to  review  the  progress  and  stimulate 
the  dialogue  between  researchers  and  users.  Copies  of  three  reports  covering 
most  of  our  recent  work  were  distributed. 

A wide  cross  section  of  people  from  industrial,  educational  and  govern- 
ment organizations  attended  the  symposium  (see  attachment  1),  reflecting  a 
widespread  interest  in' this  research  program.  Many  of  the  attendees  expressed 
desire  (1)  to  get  more  information  on  our  process  modeling  activities,  (2) 
to  get  a copy  of  our  process  modeling  computer  program  SUPREM  I,  (3)  to  be 
on  our  mailing  list  to  get  information  on  the  future  progress.  We  have  thus 
far  disseminated  28  copies  of  SUPREM  I on  magnetic  (see  attachment  2)  tape 
and  have  a dozen  additional  requests  for  copies. 

All  of  the  presentations  were  recorded  on  video  tapes  for  later  use. 
Because  of  the  immense  interest  shown  by  several  industrial,  educational, 
and  government  organizations,  copies  of  the  video  tape  will  be  shortly  avail- 
able for  wide  distribution.  This  will  be  highly  useful  for  those  who  could 
not  attend  the  symposium  and  for  internal  use  at  industrial  locations. 

During  and  after  the  symposiim,  discussions  took  place  between  the 
researchers  working  on  this  program  and  the  attendees,  which  brought  out 
interesting  points.  People  remarked  that  they  have  rarely  seen  such  an 


84 


integrated  research  program  involving  people  from  disciplines  and  organi- 
zations. Interesting  technical  questions  were  raised,  which  will  be  helpful 
to  us  in  the  future. 
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